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Abstract. We study the coupled dynamics of two populations of random replicators by means of statis- 
tical mechanics methods, and focus on the effects of relative population size, strategy correlations and 
heterogeneities in the respective co-operation pressures. To this end we generalise existing path-integral 
approaches to replicator systems with random asymmetric couplings. This technique allows one to formu- 
late an effective dynamical theory, which is exact in the thermodynamic limit and which can be solve for 
persistent order parameters in a fixed-point regime regardless of the symmetry of the interactions. The 
onset of instability can be determined self-consistently. We calculate quantities such as the diversity of the 
respective populations and their fitnesses in the stationary state, and compare results with data from a 
numerical integration of the replicator equations. 
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1 Introduction 

The evolutionary dynamics of populations of interacting 
species is often described by so-called replicator equa- 
tions (RE) [2(5]. Each species here carries a fitness, and 
species who are fitter than the average increase in con- 
centration whereas the weight of species less fit than av- 
erage decreases. RE have also found applications in game 
theory and economics [3], and an equivalence to Lotka- 
Volterra equations of population dynamics can be estab- 
lished. In the context of game theory RE describe situ- 
ations in which a game is played repeatedly by players 
chosen randomly from populations of agents. Each indi- 
vidual agent here plays one fixed (so-called pure) strategy 
throughout his lifetime, and reproduction occurs in pro- 
portion to the payoff accrued. Offspring then inherit the 
strategy of their parent. In this situation the replicators 
are hence the pure strategies of the game under consider- 
ation. In this paper we will use the interpretations of RE 
in population dynamics and in evolutionary game theory 
simultaneously, and will hence refer to the replicating en- 
tities as 'species' and 'strategies' synonymously. Similarly, 
'fitness' and 'payoff' will be used interchangeably. 

The analysis of eco-systems and replicators with fixed 
random couplings was initiated by the seminal work of [4] , 
and although the study of real-world eco-systems has since 
then moved on the incorporate more realistic and evolv- 
ing networks (see e.g. [5] and references therein), replicator 
systems with quenched random couplings are still of inter- 
est from the statistical mechanics perspective and as an 
analytically tractable benchmark. Techniques from disor- 



dered systems theory have first been applied to such mod- 
els in [HEmE] an d focus both on static approaches based 
on the replica method and on dynamical studies by means 
of generating functional techniques. The latter approach 
here has the advantage of being able to address couplings 
of general symmetry [TlllOj. whereas replica theory is re- 
stricted to cases of symmetric interaction (see [TTirnZl[T51 
[T^l[T5HTf5] for further replica studies of random replicator 
models). In the course of these studies complex non-trivial 
behaviour of random replicator systems has been identi- 
fied, with several types of phase transitions between stable 
and unstable regimes, and multiple patterns of ergodicity 
breaking. While most of the existing studies have focused 
on so-called single-population models, corresponding to 
games with only one type of player, multi-population mod- 
els have recently been addressed in [17], where a relation 
between the stability of two-population replicator systems 
and the corresponding two-player bi-matrix games has 
been considered. These are games in which the two in- 
teracting players take different roles (e.g. female versus 
male), and may hence have different sets of strategies at 
their disposal (as opposed to single-population games such 
as rock-papers-scissors, where the strategy pool is iden- 
tical for both players). Static studies of random matrix 
games can also be found in [18, 151120]. 



The aim of the present work is to extend the studies of 
[T7] to the case of two interacting populations of different 
relative size and to incorporate correlations between pay- 
off matrices as well as heterogeneity in the co-operation 
pressures of species. We first present the details of the 
model in the next section, and then formulate an effective 
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macroscopic theory in the thermodynamic limit. Specific 
results are then reported in section [U where we test the 
predictions of the resulting fixed-point theory against nu- 
merical simulations. A summary and outlook conclude the 
paper. 



2 Model 

We consider two populations of replicators (abbreviated as 
PI and P2 in the following), where we denote the species 
in system PI by i = 1, . . . , Ni, and the ones in population 
P2 by j = 1, . . . , N% . Ni and N2 are hence the sizes of the 
two populations (equivalently the number of pure strate- 
gies available to the two types of players of the bi-matrix 
game). We will consider the general case Ni ^ N2, and 
will write Ni = aN and N2 = CX2N in the following. The 
statistical mechanics analysis is then concerned with the 
limit N,Ni,N2 — > 00, where the ratios a% — Ni/N and 
a?2 = N2/N remain finite. The aspect ratio a = a.i/a.2 is 
hence a control parameter of the model. 

We will refer to the species in PI as species of type X, 
and to species in P2 as species of type Y. The composition 
of PI at time t is then characterised by a concentration 
vector (xi(i), . . . ,0^(4)), where Xi indicates the relative 
weight of species i in PI. Similarly, the configuration of P2 
can be described by (yi(t), . . . , UN 2 (t))- 111 the context of 
evolutionary game theory, vectors of this type characterise 
mixed strategies [3], with e.g. yj being proportional to the 
probability of a player of type Y using pure strategy j G 
{1, . . . , N2}. We will in the following use the normalisation 
N^ 1 J2i Xi = 1 and N^ 1 yj — 1, which will allow us 
to formulate a well-defined thermodynamic limit of the 
problem. 

Species of type X are assumed to interact only with 
species of type Y and vice versa, and the populations are 
taken to follow the following replicator equations 



V] = Vj 



-2uxXi + ^2 a ijVj ~ h 



' N 2 \ 

-2u 2 yj + ^ b Ji x i ~ h J • 

v i=l J 



(1) 



U\ and U2 are here the co-operation pressures acting on 
the two populations. Depending on their strengths u\ and 
U2 limit the growth of individual species and drive PI 
and P2 towards states with many surviving species and 
high diversity [2], We will mostly consider the situation 
in which u\,U2 > 0. The payoff matrices (or inter-species 
couplings) {o,ij,bji} are drawn from a Gaussian distribu- 
tion with the following first two moments: 



= bn = 0, 



The overbar denotes an average over the disorder, i.e. over 
samples of the random payoff matrices. All other covari- 
ances are taken to vanish. The scaling of the covariance 
matrix with the system size N is here again chosen to 
guarantee a well-defined thermodynamic limit N — > 00, 
with which the statistical mechanics theory will be con- 
cerned. P G [—1,1] is a model parameter, characterising 
the symmetry or otherwise of the interaction: for P = — 1 
we have ay = — bji corresponding to prey- predator rela- 
tions in population dynamics, and to a zero-sum game in 
the context of evolutionary game theory. For P = a,j 
and bji are uncorrelated, while P = 1 corresponds to the 
case of symmetric interaction, a^- = bji . Intermediate val- 
ues of P allow for a smooth interpolation between these 
cases, fx and /2 finally are the time-dependent mean fit- 
nesses of species in the two respective populations, 
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i=i 
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where f± = —2u\Xi + J2j a ijVj denotes the fitness of 

species i (the superscript indicates that this species be- 
(2) 

longs to PI), fj is defined similarly. These definitions 
ensure that the replicator equations (fTJ) conserve the nor- 



Ini- 



malisations N l 1 J2i2i x i — N 2 1 J2j Vj — 1 m time, 
tial conditions in our simulations are chosen to respect 
this normalisation. 



3 Statistical mechanics theory 
3.1 Uncorrelated strategies 

The dynamics of the model can be reduced to two coupled 
stochastic processes for a representative pair of 'effective' 
species, one from each population. We will not detail the 
mathematical steps here, but will only quote the final out- 
come. Similar calculations can be found in [Til 10] or, in 
a different context, in the textbook [5T]. The generating 
functional analysis delivers the following effective dynam- 
ics: 



- 2u x x{t) + ra 2 J dt'G2(t,t')x(t')+r) 1 (t)- fi{t) 



- 2u 2 y(t) + rai / dt'dit, t')y(t') + m (t) - f 2 (t) 



±{t) = x(t) 



m = y(t) 



where r)i(t) and 772 (*) represent coloured Gaussian noise 
of zero mean and with covariances to be determined self- 
consistently as 



(4) 



Hj bji — N . 



(2) 



(m(t)vi(t')) = a 2 (y(t)y(t')) s 

(V2(t)V2(t')) = a! (x(t)x(t')) 



(5) 
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(•) here refers to an average over realisations of the effec- 
tive processes, i.e. over paths of The response 
functions G\(t,t') and G 2 (t,t'), in turn, are given by 



r (t t>) / Sx{t) 



G 2 {t,t') 



) 

I ML 



(6) 



respectively. As usual in the dynamics of disordered sys- 
tems these processes are non-Markovian (as reflected by 
the retarded interaction terms), and the resulting single- 
particle noises {77! (t) , 772 (t) } exhibit non-trivial temporal 
correlations. The description in terms of the above effec- 
tive processes is equivalent to the original problem in the 
thermodynamic limit in the sense that combined disorder- 
species-averages in the microscopic problem can be eval- 
uated as averages over realisations of the effective single- 
particle noise (see 21 for further technical details regard- 
ing the generating functional technique). 

Following [7j, we proceed with a fixed-point analysis 
of ([¥]), based on the assumption of time-translation invari- 
ance and finite integrated responses 



[TO] 1 /qi and l/q2 serve as measures of the resulting di- 
versities as well, qi is here given by qi — N^ 1 J2 i xf at 
the fixed point, and a similar definition of 172 applies. If 
e.g. all species are present at equal concentration xi = 1 
in PI (i = 1, .. . ,Ni), then l/qi = 1. If however only a 
few species survive in Pi, then l/qi — > (in the extreme 
case of only one species surviving, say x\ = N\ and Xi = 
for i > 1 one has qi — Ni which tends to infinity in the 
thermodynamic limit). A detailed inspection shows that 
l/qi and l/q2 are closely related to what is referred to 
Simpson's index of diversity in ecology [22] . 



3.2 Correlated strategies 

We will below also consider the case of correlated strate- 
gies as proposed in [18,19]. Here, the moments of the ma- 
trix elements {a^ , 6^} are further constrained by impos- 
ing 

ci 1_ 

N N' 

C2 1_ 

N N ' 



Qij ^kj 



(9) 



X,. 



lim 



d,TG n (t + r, t) < 00 



(7) 



(n = 1,2). We also write q\ — (x 2 ) and q 2 = (y 2 ), 
where x and y denote the fixed point values of the single- 
species concentrations obtained from ([4]) as t — > 00. Since 
each realisation of the single-species trajectories is a time- 
dependent stochastic process, the resulting fixed point val- 
ues x and y are static random variables. See [7][T0] for 
further technical details regarding this approach. 

The fixed-point ansatz then leads to the following six 
equations for the persistent order parameters 
{<7i,92,Xi,X2,/i,/2}: 



Xi(2wi 


- a 2 rx 2 ) 


= goW, 






X2(2w 2 


- «iPxi) 


= 90(^2), 


X = X 


— 2u±x 


{a 2 q 2 )- 1 ' 2 {2u 1 


- a 2 rx 2 ) 


= 9iW, 






(a 1 q 1 )- 1 / 2 (2u 2 


- aiPxi) 


= 9i(A 2 ), 


y = y 


- 2u 2 y 


qi/(a 2 q 2 )(2ui - 


- "2PX2) 2 


= 92(Ai), 






q 2 /(axqi)(2u2 - 


- ai-^Xi) 2 


= .92(^2), 


(8) 





with correlation parameters c\,c 2 > 0. We will further- 
more restrict the analysis to the case P = — 1 when con- 
sidering strategy correlations, we then have bji = —aij. c± 
and c 2 measure the correlations within rows and columns 
of the payoff matrix. I.e. if c 2 is large, the elements 
a\j, a 2 j, a^j, . . . , tZjVij bear some correlation for any fixed 
j, reflecting a reduction of the variability of strategies in 
PI (all pure strategies i = 1, 2, 3, . . . , Nx in PI bear some 
degree of similarity). Similarly, an increased value of c\ 
makes pure strategies in P2 more alike, and hence reduces 
the strategic options of individuals in that population. 

These strategy correlations can be incorporated in the 
path-integral analysis, and lead to the following modifica- 



2u lX - / dt'G 2 (t,t')(x(t')+c 2 ) + r) 1 (t)- h(t) 



where A\ = —fi/s/ ct2Q 2 and A 2 = —fe/y/ctiqi an( i where 
g k (A) = -^= e - z ^ 2 {A - z) k for k e {0,1,2}. These 
equations are readily solved to yield the persistent order 
parameters {91,52, Xi: X2, fi,h} as functions of the model 
parameters {u%, u 2 , a±, a 2 , P}. It is here efficient to pro- 
ceed as proposed in [15] and to obtain parametric solutions 
in terms of {Z\i,Z\2}, i.e. to fix the values of these two 
quantities, and then to obtain {qi, q 2 , Xi: X2, u\, 1*2} from 
the left-hand sides of ([5]). Note that <pi = ga(Ai) and 2 = 
(70(^2) are the fractions of surviving species, in the two re- 
spective populations. In the context of game theory these 
expressions correspond to the fractions of pure strategies 
played with non-zero probabilities [18, 19, 20J. As seen in 



(10) 

The covariances of the effective single-species noises are 
now given by 



ivi(t)vi(t')) = (y(t)y(t r )) + Cl , 

(V2(t)v 2 (t')) = (x(t)x(t'))+c 2 , 



(11) 



and {ili{t)-q 2 {t')) = as before. The resulting equations for 
the persistent order parameters undergo the correspond- 
ing changes as well (we do not report them here). 



3.3 Stability 

The stability of the fixed point identified and used in the 
above ansatz can be investigated by means of a linear 
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Fig. 1. (Colour on-line) Phase diagram of the model with two 
co-operation pressures u\ and U2 (qi = 012 = 1, Ci = C2 = 0). 
The lines show the onset of instability, with the stable phase 
above the respective lines. Phase boundaries are shown for r = 
1, 0.5, 0, —0.5 from top to bottom. 



expansion about the assumed fixed point values of the 
concentrations of the effective pure strategy frequencies x 
and y and of the stationary noise variables 771 and 772 . We 
will not go into details here, as the analysis follows the 
lines of |7J. The final outcome is that the system is stable 
whenever 

—n— o > aia 2 , (12) 
XI X2 

and that non-trivial and non-decaying fluctuations may 
arise whenever this condition is violated, so that the above 
fixed-point ansatz breaks down and Eqs. (0) can no longer 
be expected to describe the behaviour of the model accu- 
rately. For ai = a 2 = 1 the above condition reduces to 
the one reported in [T7j. We also remark that, if addition- 
ally u\ = U2 = u, Eqs. (0 reduce to the one-population 
results of [TlllOj. The onset of instability then occurs at 
u c = (1 + r)/(2 v / 2), as first reported in [7]. 



4 Results 

4.1 Populations with different co-operation pressures 

We first study the case of equally large populations, a± = 
a 2 = 1, and focus on the effects of the co-operation pres- 
sures Mi and U2 on the stability and behaviour of the sys- 
tem. The resulting phase diagram has been reported in 
p"7] . and is re- iterated in Fig. Q] for completeness. The 
replicator dynamics, for any given realisation of the cou- 
pling matrices, evolve into a unique stable fixed point 
at large co-operation pressure. Below the phase transi- 
tion lines depicted in Fig. [TJ the system becomes non- 
ergodic in the sense that multiple microscopic stationary 
states exists, and initial conditions determine which of 
these is assumed asymptotically. For fully symmetric cou- 
plings r = 1 the system still evolves into a fixed point, but 
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Fig. 2. (Colour on-line) Diversities 1/qi and l/q2 of the two 
populations as a function of the co-operation pressure in pop- 
ulation 2. ui = 1 fixed (on =0.2 = l,ci = 02 = 0. Symbols are 
from simulations for TV = 100 (5000 integration steps, 10 sam- 
ples), with black circles corresponding to and red squares 
to 1/(72. r = 1, 0, —1 in the top, middle and bottom panel re- 
spectively. Lines are from theory (solid in stable phase, dashed 
in unstable regime). 



exponentially many (in N) fixed points can be expected 
to exist 6,7,8. For r < 1, i.e. for systems with partially 
or fully uncorrelated or anti-correlated couplings, no such 
behaviour is observed, and the dynamics may evolve to- 
wards a volatile, potentially chaotic state. 

The role of the two co-operation pressures is to drive 
the respective populations into the internal region of their 
concentration simpliccs (defined by the constraints 

Ni J2i=i x i = 1 an d ^2~ 1 Z~2j Vj = !)• Hence co-operation 
pressure promotes the survival of many pure strategies, or 
in other words a diverse stationary state. In order to study 
the effects of u\ and U2 we report results for a fixed value 
m = 1 in Fig. [5]and depict results for the diversity param- 
eters l/gi and 1/(72 of the respective populations, as U2 is 
varied. As expected both diversities are increasing func- 
tions of U2- The effect on population 2 is here much more 
direct, but remarkably, for small values of U2 the diversity 
of population 1 is a steep function of U2 as well, at least 
for symmetric and asymmetric couplings. This indicates a 
feedback from population 2 onto population 1: decreasing 
the co-operation pressure on population 2 leads to deple- 
tion in this population and this may then impact on the 
diversity of the other population as well, even though u\ 
is kept constant. 

In Figs. [3] and 2] we report on the fitness of the two 
populations. One here has to distinguish between the con- 
tributions of co-operation pressure and of the direct in- 
teraction to the overall fitness. While the Lagrange pa- 
rameters f\ t 2 entail the effects of co-operation pressures 
as well, we study the fitnesses v\ = N^ 1 Y^ij XiOiJyJ and 
V2 = N^ 1 ^2jiVjbjiXi purely due to direct interaction in 
Fig. HI As shown in Fig. [3] an increase of U2 typically re- 
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Fig. 3. (Colour on-line) Lagrange multipliers f\ and fa of 
the two populations as a function of the co-operation pressure 
in population 2. ui = 1 fixed («i = oli = l,cx = C2 = 0). 
Simulation parameters, lines and symbols as in Fig. [2] Circles 
are population 1, squares population 2. 



duces f\ and f 2 , provided P > 0, i.e. provided that there 
is no degree of anti-correlation present in the coupling ma- 
trices. Again the effects of a change in u 2 are much more 
pronounced in P2 than in PI, and f± remains almost con- 
stant except for a region at small values of u-i- For negative 
values of P the effect is quite different, as seen in the lower 
panel of Fig. [3] While, as before, is a decreasing func- 
tion of U2 , non-trivial behaviour is induced in population 
1, and /i is found to be increasing in u%, especially when 
U2 small. This is the regime in which ui strongly controls 
the diversity of population 2, and a reduction of the di- 
versity in P2 appears to impact negatively on the fitness 
of PI. 

In order to disentangle the effects of co-operation pres- 
sure related contributions to the fitnesses, we study the 
behaviour of v n = f n + 2u n q n , n = 1,2 in Fig. 0J For 
positive strategy correlation, we find that both fitnesses 
decrease as 112 in increased. The non-trivial behaviour at 
negative strategy correlation remains, however, and v\ is 
found to be a non-monotonous function of U2, again sig- 
nalling an indirect impact of the diversity of P2 on the 
payoff of PI. 

An interesting effect is observed for the case P = — 1, 
i.e. for full anti-correlation, see Fig. [5J In order to char- 
acterise the behaviour of 1^2, we have here extended the 
range of the co-operation pressures to include small neg- 
ative values of U2- In the zero-sum case P = — 1 one has 
v\ = —V2 by construction, and both fitnesses are found 
to be non-monotonic functions of 1*2 at fixed co-operation 
pressure in population 1. In particular a maximum of the 
fitness of P2 is found as the co-operation pressure U2 is in- 
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Fig. 4. (Colour on-line) Mean fitness V\ and V2 of the two 
populations as a function of the co-operation pressure in pop- 
ulation 2. ui = 1 fixed (ai = 02 = 1, ci = C2 = 0). Simulation 
parameters, lines and symbols as in Fig. [2] but P = —0.5, 0.5 
are chosen here. Symbols are from simulations, N = 100 (5000 
iteration steps, 20 samples) with circles marking the fitness of 
population 1, squares the one of population 2. 



creased, whereas the other population experiences a min- 
imum in fitness under these conditions Q 



4.2 Effects of relative population size 

The effects of varying the relative system size of the two 
populations are illustrated in Figs. [6] and [7] We here first 
follow 18,19 and fix the size parameter a\ in the first 
population, and vary the size of P2. As shown in Fig.[6l the 
diversities of both populations decrease as a<2 = N%/N is 
increased. The effect is much stronger in PI than in P2. In 
Fig. [7] we depict the fitnesses Ui, V2 of the the two popula- 
tions for this case of rectangular payoff matrices a\ 7^ ot2- 
We find that, within our settings and for P > 0, both 
populations profit from reducing the aspect ratio ax/a^, 
where the increase of fitness appears more pronounced 
for PI than for P2. Thus it appears it is the size of the 
population with which a given species interacts, rather 
than the size of its own population, which determines 
the fitness of this given species. Inspection of the case 
of full anti-correlation (lower panel in Fig. [7J shows that 
the fitness of PI may display non-monotonic behaviour 
as a function of ot2 ■ We note here that v\jv% = 0,2 for 
the case of fully symmetric couplings (P = 1), and that 



1 We here remark that while numerical experiments agree 
near perfectly with the theoretical predictions for positive U2, 
they become less reliable at U2 ~ 0, where the theoretical 
curves are fairly steep. The accuracy of simulations might here 
suffer from a small number of surviving species and the asso- 
ciated finite-size effects and sample-to-sample fluctuations. 
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Fig. 5. (Colour on-line) Mean fitness and z/2 of the two 
populations as a function of the co-operation pressure in pop- 
ulation 2. Mi = 1 fixed. Simulation parameters, lines and sym- 
bols as in Fig. [2] but r = — 1. Symbols are from simulations, 
TV = 100 (5000 iteration steps, 200 samples) with circles mark- 
ing the fitness of population 1, squares the one of population 
2. 

v\jv2 = —ot2 for F = — 1, i.e. for a%j = —bji. These re- 
lations are due to the construction of the model, where 
one has v\ = TVf 1 £V x l a lj y ] and v 2 = N^ 1 J2ji UjbjiZi 
as well as N1/N2 = a 1/0:2. No such simple relations be- 
tween the resulting v\ , v 2 are found for — 1 < F < 1 , as 
dij and bji are then neither fully correlated nor fully anti- 
correlated. 

The results described so far in this section are not 
found to be symmetrical around 02 = 1, In a symmetri- 
cal situation, where only the aspect ratio 01/02 matters, 
one would expect that choosing 02 = a is equivalent to 
setting 02 = a -1 for any a > 0, subject to a relabelling 
of the two populations. Recall, however that the coupling 
strengths afj — b 2 ^ remain fixed and are given by l/N. 
The respective sizes of the two populations are Ni = a\N 
and N2 — 02-/V, where we fixed Ni = N above. These def- 
initions follow those of [18pi9j . A more symmetrical setup 
can be achieved by setting e.g. Oi = y/a and 02 = 1/y/a, 
with < a < 1 the aspect ratio a = 01/02 = Ni/N%. 
As seen in Fig. the diversity of population 1 decreases 
as its relative size o is increases, whereas I/92 is mono- 
tonically decreasing (we do not report simulation results 
here to keep the figure readable). The behaviour of the 
fitnesses crucially depends on the symmetry of couplings, 
with multiple crossings of V\ and v 2 observed for suitable 
model parameters as shown in Fig. [8] 



4.3 Strategy correlations 

Finally, we have examined the effects of strategy corre- 
lations. We here fix oi = 02 = 1 and F = — 1. Results 
are presented in Fig. [5J where we depict the case of vary- 
ing C2 at fixed c\ — 0. As seen in the figure, the diversity 
(or equivalently the fraction of pure strategies played with 
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Fig. 6. (Colour on-line). Diversity parameters 1/qi (circles) 
and l/q2 (squares) as the size «2 of P2 is varied (on = 1 
remains fixed; Ui = u% = 1, ci = c% = 0). Simulations are for 
N — 200, 20000 integration steps, 10 samples, solid lines are 
from the theory. F — 1,0, —1 from top to bottom. 



non-zero probability) decreases in both populations as C2 
is increased. The effect is stronger in P2 than in PI. At the 
same time, the payoff for population two increases as the 
correlation parameter C2 is increased, whereas the payoff 
for PI decreases (note that v\= —v^ by construction, as 
we are considering the case of zero-sum games, P = —1). 
The condensation of P2 into a small subset of strategies 
can here be understood as follows: an increased value of 
C2 induces correlations within the columns of the payoff 
matrix A = (ay), so that in the extreme case of large 
C2 S> I, the {ay}j=i,...,iv become mostly independent of 
i. Thus the payoff for both players becomes more and more 
independent of the action of player X, and depends only 
on the choice j of player Y. In other words, the strategies 
available to player X become more and more alike as C2 is 
increased, thus rendering some of Y's strategies generally 
beneficial for Y (independently of the choice of X), and 
others detrimental. Hence, Y will mostly play strategies 
beneficial from his perspective and will reduce the diver- 
sity 1/92 of his actions, while at the same time increasing 
his payoff as seen in Fig. [9l See also [18] for similar 
cases. A more symmetrical situation can be constructed 
by considering c\ — C2 = c, a case in which correlations are 
present both in the rows and columns of the payoff matrix. 
We do not depict results here, but only remark that in this 
case the diversity of both populations decreases as c is in- 
creased (the payoff for both populations remains strictly 
zero in this case due to the zero-sum property of the games 
considered, and the exchange symmetry between PI and 
P2). 
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Fig. 7. (Colour on-line). Fitnesses v\ (circles) and V2 (squares) 
of the two populations, as the size ai of P2 is varied (ai = 1 
remains fixed; m = 112 = 1). Simulations are for TV = 200, 
20000 integration steps, 10 samples, solid lines are from the 
theory, P — 1, 0, —1 from top to bottom. 
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Fig. 9. (Colour on-line). Effects of strategy correlations. c\ is 
fixed to ci = as C2 is varied from to 1. Fully anti-symmetric 
(zero-sum) case F = — 1. Co-operation pressures are ui = «2 = 
1. Results from the fixed-point ansatz are shown as solid lines, 
markers are from simulations, with circles corresponding to 
PI and squares to P2 (N = 100, 10000 integration steps, 200 
samples) . 




Fig. 8. (Colour on-line). Diversities and fitnesses v\ of the two 
populations as a function of the aspect ratio a (qi = 

Upper panels are for P = —0.5, lower panels for r = 



-1/2-1 



(ui 



«2 



1). Solid lines are PI, dot-dashed lines P2. 



4.4 One-population models with heterogeneous 
co-operation pressure 

In this section we will consider a single-population replica- 
tor model with heterogeneous (i.e. species-dependent) co- 
operation pressures. Specifically consider N species sub- 
ject to the replicator equations 



N 

'ln,.r, + w ij x o f I 



.,N (13) 



where the coupling matrix is a Gaussian random quantity 
as in [7] (Wpj — 0,w^ = l/N,Wi]w]i = r/N), but where 
the co-operation pressure now carries a species index i, 
and where each itj is assumed to be drawn independently 
from the distribution p(u) at the beginning of the dynam- 
ics and then remains fixed. Note that while the overall 
concentration iV -1 J2 i Xi = 1 is conserved, species subject 
to a certain (e.g. high) co-operation pressure may well die 
out or be reduced in concentration to the advantage of 
species of a lesser co-operation pressure. 

The further analysis leads to an ensemble of effective 
processes, one for each value of u in the support of p(u): 



■2„x n P I dt'G{t 1 t l )x u {t')+r lu {t) - f(t) . 

(14) 



The response function is now defined as 
G(t,f) = J dup{u) (Sx u (t)/Srj u (t')). Furthermore we have 
(i]u(t)r]u(t')} = J du'p(u') {x u i{t)x u i (t')) independently of 
u. A fixed point ansatz results in the following self-consistent 
equations for the integrated response x> the (inverse) di- 
versity parameter q and the fitness /: 



X 

-1/2 
1 



9o(A) 
92(A) 



du 
du 
du 



p(u) 



(2u - r x ) ' 

p( u ) 
(2u-r x y 

p( u ) 
\2u-r x ) 2 ' 



(15) 



where A = —f/^jq. While we note that the fitnesses / of 
any surviving species come out as equal (and independent 



8 



Tobias Galla: Two-population evolutionary dynamics 



of their co-operation pressures), their relative concentra- 
tions 

(16) 



c[u) = W)\ ^ 



and second moments 



Q(u) = 



I(u)\ 



(17) 



(suitable sample-averages are implied) maybe well be de- 
pendent on u Q. Here I(u) denotes the set of species i e 
{1, . . . , N} with co-operation pressure u t € [u— du, u+du], 
with du an infinitesimal element. Overall self-consistency 
requires that J dup(u)C(u) — 1 and j dup{u)Q(u) = q. 

We will here restrict to the case of a flat distribution 
p(u) over the interval [p — s, p + s], so that the integrals 
on the right-hand side of (|15|) can be performed explicitly 
to give 



X = 9o(A)(A S )-Hn 

=g 1 (A)(4s)- 1 \n 
92(A) 



2p + 2s- F X 
2p-2s- r x 
2p + 2s - r x 



2p-2s- r x 



[2p - r x ) 2 - 4s 2 



(18) 



Results are shown and compared with simulations in Fig. 
[TUl Interestingly one observes a decline in diversity of the 
population as the variability s of the co-operation pressure 
is increased. Thus increasing the complexity of the prop- 
erties of the individual species might lead to a less diverse 
composition of the eco-system at stationarity. This effect 
seems to be mostly independent of the symmetry of the 
coupling^. The inset of Fig. \W\ demonstrates that the 
relative weight C(u) of species subject to co-operation u 
decreases non-linearly with u, so that species with compa- 
rably low co-operation pressure dominate the population 
at the fixed point. More precisely C(u) is found to be 
given by C(u) = y/qgi (A)/(2u — r X ), i.e. it roughly de- 
cays as the inverse power of u. As seen in the inset of the 
figure the theoretical prediction of this behaviour agrees 
perfectly with results from simulations. 



5 Conclusions 

In summary, we have extended existing generating func- 
tional techniques to study the behaviour of two-population 
replicator systems, and have focused on the effects of co- 
operation pressure, relative population size and strategy 
correlation. A phase transition between a stable phase 



2 We here note that the fraction of survivors 4>{u) — 
Sig/( u ) (with &(■) the step function) comes out as in- 
dependent of u, and is given by <j>(u) = tf> = go(A). 

3 We have not been able to solve the resulting fixed-point 
equations for _T = 1 at large widths s of the distribution of 
co-operation pressures. 




Fig. 10. (Colour on-line) One-population model with hetero- 
geneous co-operation pressures. Main panel: diversity versus 
variability of co-operation pressure (p = 1 fixed). Upper curve 
shows r = — 1, lower curve f = 0. Symbols from simulations 
(N = 300, 5000 integration steps, 10 samples). Inset: concen- 
tration of sub-population with co-operation pressure u as a 
function of u (r = 0, p = 1, s = 0.5). Solid line is theory, noisy 
line simulations (N — 300, 100 samples). Dashed line to the 
lower right in main panel indicates unstable regime for r — 



with one unique fixed point of the replicator dynamics 
has been identified, and characterised analytically. This 
phase is separated from a second unstable phase by a tran- 
sition line in parameter space, which can be determined 
from the statistical mechanics analysis. Our study demon- 
strates that control parameters such as co-operation pres- 
sure, aspect ratio and strategy correlation have non-trivial 
effects on the the system of replicators, and can induce 
non-monotonic behaviour of the resulting fitnesses. We 
have also addressed single population models with species- 
dependent co-operation pressure. The statistical mechan- 
ics theory is then formulated in terms of an ensemble of 
effective processes, one for each value of the co-operation 
pressure present in the population. We find that variabil- 
ity in the co-operation pressures reduces the diversity of 
the set of surviving species. In such a mixed population, 
the weight of species subject to a specific co-operation 
pressure u scales as w -1 asymptotically. 

Natural extensions of the present model include the 
generalisation to a larger number of populations, and the 
study of dynamics which is different from standard repli- 
cator equations [Tll^l[^5]. Furthermore, relatively little is 
known about the non-ergodic phase of random replicator 
systems, so that future work might address the proper- 
ties of such phases. From the point of view of real- world 
eco-systems and population dynamics the assumption of 
a fully connected interaction matrix is at best a crude 
approximation, and would ideally need to be replaced by 
ensembles with sparse interactions. While the analysis of 
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fully connected replicator models like the one discussed 
in this paper, is relatively straightforward and leads to 
an effective theory in terms of two-time quantities (the 
correlation and response functions), dilute replicator sys- 
tems pose a much more demanding challenge, as order 
parameter equations do not close on the two-time level 
[24 . Future work might hence address such models, po- 
tentially relying on recently developed techniques to study 
spin systems with sparse interaction matrices |25j . 
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